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ABSTRACT 

A  computational  chemistry  toolkit  is  presented, 
demonstrating  the  ability  to  predict  properties  related  to 
performance  or  hazard  of  these  materials.  The  purpose 
of  this  toolkit  is  to  allow  the  design  and  assessment  of 
advanced  energetic  materials  before  investing  in 
synthesis.  Newly  developed  methods  are  presented  to 
treat  emerging  exotic  high-nitrogen  systems  considered  as 
candidates  for  insensitive  and  environmentally-friendly 
advanced  energetic  materials. 

1.  INTRODUCTION 

New  computer  tools  and  software  developments  have 
made  possible  “Materials  By  Design”  programs,  in  which 
new  materials  can  be  computationally  constructed  and 
assessed  before  actual  production.  National  and 
international  investigators  from  almost  every  university, 
institution,  research  center,  and  within  industry  are 
aggressively  engaged  in  efforts  to  develop  multiscale 
modeling  and  simulation  methodologies  to  model  and 
analyze  complex  phenomena  across  multiple  time  and 
length  scales.  It  is  now  possible  to  predict 
structure/property  relationships  for  both  simple  and 
complex  materials,  as  well  as  predict  the  dynamic 
behavior  of  materials  at  all  relevant  length  scales, 
including  the  more  challenging  and  least  empirical 
atomistic  level.  However,  only  in  recent  years  have 
fundamental  theoretical  chemistry  modeling  and 
simulation  tools  for  energetic  materials  (EM)  research 
dramatically  evolved.  These  efforts  are  expected  to 
impact  the  DOD’s  ability  to  model,  analyze,  and 
understand  a  vast  array  of  energetic  materials  and 
complex  formulations,  leading  to  a  Materials  By  Design 
capability  in  this  area. 

It  is  not  sufficient,  however,  simply  to  collect 
computational  models  that  predict  various  properties  and 
processes  of  interest  in  a  notional  material.  Rather,  a 
protocol  and  framework  for  integration  and  validation  of 
emerging  methods  and  models  for  energetic  materials 
design  is  required.  We  have  designed  such  a  framework 
within  a  Computational  Chemistry  Toolkit  for  energetic 
materials  design.  It  is  hoped  that  this  toolkit  will  evolve 
into  a  collection  of  well-integrated  multiscale  modeling 
methodologies  capable  of  coupling  behaviors  from  the 
atomic  scale  through  to  the  frill-scale  system.  At  this 
time,  the  majority  of  our  efforts  have  been  devoted  to 


atomic-level  modeling,  which  will  be  the  focus  of  this 
paper. 

Significant  advances  have  been  made  in  the 
development  of  computational  atomistic  methods  to 
predict  properties  related  to  performance,  sensitivity,  and 
environmental  impact  of  a  new  material.  Special 
attention  has  been  given  to  the  prediction  of  two 
properties  that  are  used  to  provide  an  initial  assessment  of 
the  potential  performance  of  a  material  in  a  gun  or 
warhead:  the  heat  of  formation  and  the  density  of  the 
material.  Performance  metrics  (impetus  or  detonation 
velocities  and  pressures)  are  dependent  on  the  energy 
content  of  the  charge,  reflected  by  the  heat  of  formation 
of  the  energetic  material,  and  the  density,  which  is  an 
indicator  of  how  much  material  can  be  packed  into  the 
charge.  Therefore,  a  high-density  material  with  a  high 
heat  of  formation  would  appear  to  be  a  desirable 
candidate.  We  have  tested  and  assessed  a  variety  of 
computational  chemistry  methods  to  accurately  predict 
such  performance  properties  for  conventional  CHNO 
explosives  (Rice,  Pai,  and  Hare,  1999;  Byrd  and  Rice, 
2006;  Rice  and  Hare,  2002;  Rice  and  Sorescu,  2004). 
Unfortunately,  until  recently  there  has  been  no  similar 
collection  of  computational  methods  available  to  predict  a 
priori  these  properties  for  high  nitrogen  compounds  or  for 
ionic  crystals  (either  high-nitrogen  or  CHNO  salts). 
Previous  tools  that  have  been  developed  for  CHNO 
containing  energetics  are  completely  unsuited  for  ionic 
materials,  and  cannot  be  assumed  to  be  transferable  to 
high-nitrogen  neutrals.  Therefore,  in  such  cases,  new 
predictive  methods  are  required.  By  expanding  on  earlier 
methodologies  (Baboul  et  al.,  1999),  as  well  as  testing 
alternative  methods  that  are  applicable  to  ionic  materials 
(Jenkins  at  ah,  2002),  we  are  now  able  to  predict  the  heats 
of  formation  and  densities  of  both  neutral  and  ionic  high 
nitrogen  materials. 

Other  aspects  of  materials  that  are  receiving 
increased  attention  are  the  negative  effects  of  munitions 
on  the  environment  and  the  cost  associated  with  both  test 
facility  clean-up  and  disposal  of  synthesis  waste.  These 
have  recently  been  introduced  into  the  factors  that  are 
considered  before  investing  in  synthesis  and  testing. 
Accordingly,  we  are  expanding  our  program  to  develop 
models  to  address  the  environmental  hazards  of  energetic 
materials.  Through  the  use  of  quantitative  structure 
property/activity  relationships  (QSPR  /  QSAR),  we  have 
begun  investigating  possible  indicators  of  mammalian  and 
aquatic  toxicity. 


1 


Report  Documentation  Page 

Form  Approved 

OMB  No.  0704-0188 

Public  reporting  burden  for  the  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instructions,  searching  existing  data  sources,  gathering  and 
maintaining  the  data  needed,  and  completing  and  reviewing  the  collection  of  information.  Send  comments  regarding  this  burden  estimate  or  any  other  aspect  of  this  collection  of  information, 
including  suggestions  for  reducing  this  burden,  to  Washington  Headquarters  Services,  Directorate  for  Information  Operations  and  Reports,  1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington 

VA  22202-4302.  Respondents  should  be  aware  that  notwithstanding  any  other  provision  of  law,  no  person  shall  be  subject  to  a  penalty  for  failing  to  comply  with  a  collection  of  information  if  it 
does  not  display  a  currently  valid  OMB  control  number. 

1 .  REPORT  DATE  2.  REPORT  TYPE 

01  NOV  2006  N/A 

3.  DATES  COVERED 

4.  TITLE  AND  SUBTITLE 

Computational  Chemistry  Toolkit  For  Energetic  Materials  Design 

5a.  CONTRACT  NUMBER 

5b.  GRANT  NUMBER 

5c.  PROGRAM  ELEMENT  NUMBER 

6.  AUTHOR(S) 

5d.  PROJECT  NUMBER 

5e.  TASK  NUMBER 

5f.  WORK  UNIT  NUMBER 

7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS(ES) 

U.  S.  Army  Research  Laboratory,  AMSRD-ARL-WM-BD  Aberdeen 
Proving  Ground,  MD  21005-5069 

8.  PERFORMING  ORGANIZATION 

REPORT  NUMBER 

9.  SPONSORING/MONITORING  AGENCY  NAME(S)  AND  ADDRESS(ES) 

10.  SPONSOR/MONITOR'S  ACRONYM(S) 

11.  SPONSOR/MONITOR'S  REPORT 
NUMBER(S) 

12.  DISTRIBUTION/AVAILABILITY  STATEMENT 

Approved  for  public  release,  distribution  unlimited 

13.  SUPPLEMENTARY  NOTES 

See  also  ADM002075.,  The  original  document  contains  color  images. 

14.  ABSTRACT 

15.  SUBJECT  TERMS 

16.  SECURITY  CLASSIFICATION  OF:  17.  LIMITATION  OF 

18.  NUMBER  19a.  NAME  OF 

a.  REPORT  b.  ABSTRACT  c.  THIS  PAGE  |J|J 

unclassified  unclassified  unclassified 

6 

Standard  Form  298  (Rev.  8-98) 

Prescribed  by  ANSI  Std  Z39-18 


This  paper  will  describe  our  past  and  current  efforts 
to  further  predictive  capabilities  of  energetic  materials. 
This  paper  will  outline  the  procedures  we  have 
investigated  in  order  to  predict  the  heats  of  formation  and 
densities  of  both  neutral  and  ionic  high  nitrogen  materials. 
We  will  also  provide  a  discussion  on  limitations  on 
current  theoretical  tools  and  subsequent  needs  and 
intentions  for  advancing  our  models. 

2.  RESULTS  AND  DISCUSSION 
2.1.  Prediction  of  Crystal  Density 

Several  approaches  have  been  used  to  predict  crystal 
densities  without  a  priori  knowledge  of  the  system.  The 
first,  and  most  sophisticated,  is  that  of  ab  initio  crystal 
prediction,  in  which  the  crystal  microstructure  is 
determined  using  information  about  a  single,  isolated 
molecule  (and  requires  a  description  of  interatomic 
interactions).  This  method  and  a  general  CHNO 
interaction  potential  were  assessed  as  viable  tools  for 
crystalline  microstructure  (and  density)  predictions  of 
explosives  by  Rice  and  Sorescu  (Rice  and  Sorescu,  2004) 
through  application  to  a  large  number  (174)  of  CHNO 
crystals  with  energetic  functionalities.  Eight-five  percent 
of  the  predictions  produced  crystallographic  parameters 
and  molecular  configurations  which  matched  those  of  the 
experimental  counterpart  and  were,  on  average,  higher 
than  experiment  by  less  than  3%.  The  overall  good 
performance  of  the  model  demonstrates  it  is  reasonable  to 
use  in  crystal  predictions  of  similar  CHNO  systems  where 
the  crystallographic  symmetry  is  not  known  a  priori. 
Unfortunately,  the  method  and  model  have  exhibited 
limitations  when  applied  to  emerging  novel  energetic 
materials. 

Using  experimental  infonnation  generated  by  T. 
Klapotke  (U.  Munich),  we  established  that  the 
aforementioned  procedure  was  unsuitable  for  prediction 
of  high-nitrogen  neutral  and  ionic  crystals.  The  main 
limitation  of  the  method  and  model  for  treatment  of  high- 
nitrogen  neutrals  was  in  the  force-field.  The  model  used 
in  the  earlier  study  was  parameterized  to  a  nitramine 
crystal,  and  shown  to  be  transferable  to  other  CHNO 
systems.  Its  application  to  high-nitrogen  systems 
provided  disappointing  predictions  of  density.  The 
method  as  currently  formulated  accepts  only  single 
molecules  as  the  initial  packing  moiety  around  which 
candidate  crystal  are  constructed.  Since  ionic  systems  are 
inherently  multi-molecular  moieties,  this  method  cannot 
be  currently  applied  to  ionics. 

In  light  of  these  complications,  we  explored  a  new 
procedure  to  predict  crystal  densities  using  quantum 
mechanical  predictions  of  the  volume  of  either  isolated 


molecules  or  ions.  In  this  procedure  the  “molecular 
volume”  (or  “formula  unit  volume”,  for  ionic  systems)  is 
used  to  calculate  the  crystal  density.  The  molecular 
volume  is  that  contained  within  specified  isosurfaces  of 
the  electron  density  of  the  molecule,  which  can  be  readily 
calculated  using  quantum  mechanics.  To  assess  this 
method  for  neutrals,  we  have  predicted  crystal  densities 
for  181  CHNO  species  for  which  experimental 
crystallographic  information  exists.  Predicted  molecular 
volumes  have  average  and  root  mean  square  (rms) 
deviation  from  experiment  -0.9  and  3.7%,  respectively. 
We  note  that  this  method  of  crystal  density  prediction 
produced  statistical  agreement  with  experiment  as  good  as 
those  generated  in  the  study  of  174  CHNO  systems  using 
ab  initio  crystal  prediction,  at  a  substantially  reduced 
computational  cost  (Rice  and  Sorescu,  2004).  We  then 
calculated  molecular  volumes  for  eight  high-nitrogen 
neutrals;  predicted  densities  have  a  2.7%  rms  deviation 
from  experiment. 

We  next  used  the  procedure  to  predict  crystal 
densities  of  ionic  crystals.  As  it  is  not  possible  to  use  gas 
phase  quantum  mechanical  predictions  of  isolated  ionic 
partners  to  determine  their  most  probable  relative 
arrangement  in  the  crystal,  (and  thus,  the  volume  of  the 
pair)  the  volume  of  the  formula  unit  (where  M 

denotes  the  cation  and  X  denotes  the  anion)  is  defined  as: 

Volume  =p\+  +  qW.  (1) 

Our  assessment  of  this  method  was  attempted  for  34  high- 
nitrogen  ionic  salts  provided  by  Klapotke  (U  Munich);  a 
comparison  of  the  predicted  and  experimental  molecular 
volumes  are  given  in  Fig.  1.  Predicted  crystal  densities 
using  this  method  had  a  4.2%  rms  deviation  from 
experiment. 


Compound  number 

Figure  1.  Experimental  and  theoretical  molecular 
volumes  for  34  high  nitrogen  salts. 

2.2.  Prediction  of  Heats  of  Formation 
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Most  energetic  materials  exist  in  the  solid  phase,  and 
their  heats  of  formation  in  this  phase  are  a  property  that 
cannot  be  directly  calculated  using  quantum  mechanical 
theory.  However  condensed  phase  heats  of  formation  can 
be  obtained  if  cohesive  energies  can  be  determined  and 
combined  with  quantum  mechanical  predictions  of  gas 
phase  heats  of  formation.  For  neutral  molecules,  the 
cohesive  energy  of  a  solid  is  roughly  that  of  the  heat  of 
sublimation;  for  ionic  crystals,  the  cohesive  energy  is  the 
lattice  energy. 

Using  a  procedure  developed  by  Politzer  (Politzer, 
1994),  we  have  parameterized  equations  that  are  used  to 
predict  the  heats  of  sublimation  and  vaporization  using 
QM  information  for  an  isolated  molecule  of  a  neutral 
crystal  (Rice,  Pai,  and  Hare,  1999;  Byrd  and  Rice,  2006) 
We  combined  this  information  with  predicted  gas-phase 
heats  of  formation  of  a  series  of  CHNO  systems  using 
atom  equivalents  (Rice,  Pai,  and  Hare,  1999).  While  this 
work  demonstrated  the  utility  of  the  method,  the  error 
associated  with  the  solid  phase  heats  of  formation 
compelled  us  to  improve  its  accuracy.  In  a  subsequent 
study,  we  expanded  our  calculations  to  use  a  larger  basis 
set  and  group  equivalents  (Byrd  and  Rice,  2006).  The 
resulting  rms  and  maximum  deviation  of  predicted  solid 
heats  of  formation  were  reduced  by  approximately  40%. 

It  could  not  be  assumed  that  these  computational 
tools,  parameterized  using  experimental  information  for 
CHNO  systems,  would  be  transferable  to  high-nitrogen 
systems.  In  fact,  application  of  the  tools  to  a  newly- 
synthesized  high-nitrogen  compound  (l-Methyl-5- 
(methylnitramino)-l//-tetrazole)  (Karaghiosoff  et  al., 
2006)  predicted  a  solid  phase  heat  of  formation  of  70.1 
kcal/mol.  far  larger  than  the  measured  value  (2.8 
kcal/mol).  Since  this  was  the  first  application  of  our  tool 
to  a  high-nitrogen  compound,  we  did  not  have  reason  to 
have  confidence  in  the  computational  tool  for  use  in 
predicting  solid  phase  heats  of  formation  for  high  nitrogen 
neutrals.  Therefore,  we  were  forced  to  try  to  identify  and 
remedy  the  most  probably  source  of  error  to  allow 
application  to  high-nitrogen  systems. 

Values  for  heats  of  sublimation  of  the  CHNO  crystals 
we  surveyed  typically  range  between  10-30  kcal/mol, 
with  very  few  systems  falling  outside  this  range. 
Additionally,  we  know  that  the  heat  of  sublimation  must 
always  be  positive.  Conversely,  experimental  gas  phase 
heats  of  formation  for  the  CHNO  systems  we  have 
surveyed  fall  within  a  much  larger  range  (from  -70  to 
+100  kcal/mol)  (Byrd  and  Rice,  2006).  In  considering  the 
disagreement  of  our  tool  with  the  measured  results  for  1- 
Mcthyl-5-(methylnitramino)-l  //-tctrazole,  we  concluded 
that  the  major  source  of  any  discrepancy  in  our 
calculations  and  the  measured  value  would  be  in  the  value 
of  the  gas  phase  heat  of  formation.  Thus,  we  focused  on 
exploring  whether  the  error  in  the  calculation  of  the  solid 


phase  heat  of  formation  might  be  due  to  improper 
prediction  of  the  gas  phase  heats  of  formation.  Therefore, 
we  pursued  other  quantum-mechanically  based  methods 
for  predictions  of  this  thermodynamic  value.  In  all 
calculations  reported  hereafter,  the  heats  of  sublimation 
were  calculated  using  the  aforementioned  method  (Byrd 
and  Rice,  2006). 

The  alternative  methods  we  are  using  for  predicting 
gas  phase  heats  of  formation  are  variants  of  the  popular 
G3  methods  specifically  the  G3(MP2)  (Curtiss  et  al., 
1999)  and  G3(B3LYP)  (Baboul  et  al.,  1999)  approaches. 
The  G3(MP2),  G3(B3LYP)  and  previously-discussed 
atom/group  equivalent  methods  were  used  to  calculate  the 
gas  phase  heats  of  formation  of  a  set  of  neutral  high 
nitrogen  crystals  and  added  to  the  predicted  heats  of 
sublimation.  The  results  are  compared  with  experiment  in 
Figure  2.  For  all  compounds,  the  theoretical  methods  are 
in  reasonable  agreement  with  each  other;  additionally,  all 
calculated  values  are  in  agreement  with  experiment 
except  for  Compound  2  (l-Methyl-5-(methylnitramino)- 
1  //-tctrazole)  (Karaghiosoff  et  al.,  2006).  This  suggests 
that  either  re -measurement  of  the  heat  of  formation  of 
Compound  2  should  be  performed  or  that  there  is 
something  within  this  system  that  is  not  captured  by  the 
quantum  mechanical  calculations.  We  emphasize  that  the 
two  quantum  mechanical  approaches  (atom/group 
equivalents  versus  G3  methods)  are  sufficiently  dissimilar 
that  it  is  extremely  unlikely  that  the  error  is  within  the 
predicted  gas  phase  heat  of  formation.  There  is  a 
possibility,  of  course,  that  the  heat  of  sublimation  is 
poorly  predicted  for  this  system.  That  would  not, 
however,  explain  the  large  discrepancy  between  the 
theoretical  values  and  the  measured  result  for  Compound 
2. 


Compound  number 

Figure  2.  Solid  phase  heats  of  formation  for  high-nitrogen 
crystals.  Experimental  values  are  taken  from 
Karaghiosoff  et  al.,  2006  and  Huynh  et  al.,  2005.  “Old 
Atom  Equiv”  denotes  calculations  using  the  atom 
equivalent  method  described  in  Rice  et  al.,  1999.  “New 
Group  Equiv”  denotes  calculations  using  the  group 
equivalent  method  described  in  Byrd  and  Rice,  2006. 
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The  calculation  of  heats  of  formation  for  ionic  solids 
is  substantially  different  than  that  discussed  previously  for 
neutral  crystals.  The  heat  of  formation  of  ionic 
compounds  is  determined  using  Born-Haber  Cycles  in 
which  a  series  of  reactions  involving  the  ionic 
components  are  employed  to  produce  the  overall  final 
result.  The  user  is  free  to  choose  any  series  of  reactions 
in  the  Born-Haber  Cycle  that  will  lead  to  the  formation  of 
the  ionic  crystal;  we  have  chosen  the  following: 

M(g)  +  X(g)  — *  M  (g)  +  X  (g)  — >  MX  (soiid)  (2) 

(AHg)  +/-  e  (AHlat) 

This  cycle  consists  of  three  steps  in  which  formation 
energies  for  each  of  the  components  can  be  calculated. 
The  first  step  (Step  1)  requires  the  evaluation  of  the  gas 
phase  heats  of  formation  of  the  neutral  form  of  the 
charged  components.  The  formation  energies  for  the 
ionic  components  are  then  determined  by  evaluating  the 
electron  affinity  or  ionization  energies  and  adding  these  to 
the  heats  of  formation  of  the  neutral  moieties  (Step  2). 
The  third  step  of  the  cycle  requires  the  evaluation  of  the 
lattice  enthalpy  (Step  3).  The  main  drawback  to  this 
approach  is  the  requirement  that  the  neutral  form  of  the 
charged  species  must  be  a  minimum  on  the  potential 
energy  surface  (PES).  If  it  is  not,  then  our  Born-Haber 
Cycle  would  require  inclusion  of  reactions  leading  to  the 
neutral  moieties  in  Step  1 . 

Unfortunately,  there  are  limited  (if  any)  experimental 
data  that  can  be  used  to  assess  the  quality  of  the 
calculations  associated  with  each  step  in  the 
aforementioned  cycle.  We  are  limited  to  comparing  the 
overall  final  result  with  the  experimental  value  of  the  heat 
of  formation  of  the  ionic  crystal.  This  makes  establishing 
the  sources  of  errors  in  this  cycle  complicated. 

For  Step  1,  a  G3-type  approach  as  described  above 
was  used  to  predict  the  gas  phase  Heat  of  Formation  for 
the  ionic  species.  Step  2,  the  evaluation  of  the  Ionization 
Energy/Electron  Affinity,  was  performed  by  calculating 
the  energy  difference  between  the  neutral  form  of  the 
ionic  component  and  that  of  the  ionized  versions.  In  both 
cases,  it  is  imperative  that  the  geometries  are  optimized, 
since  in  many  cases,  the  molecular  structures  for  the 
neutral  forms  of  the  ions  are  substantially  different  than 
the  ionic  forms.  In  evaluating  20  cations  and  6  anions 
used  in  high-nitrogen  ionic  crystals  provided  by  Klapotke 
(U  Munich)  we  found  that  half  of  the  compounds  had 
structural  relaxation  energies  on  the  order  of  50  kcal/mol. 

The  final  step  in  the  evaluation  of  the  Born-Haber 
cycle  is  the  determination  of  the  lattice  enthalpy,  a 
measure  of  the  energy  required  to  dissociate  the  ionic 
crystal  into  its  gaseous  ions.  The  magnitude  of  this 
energy  is  quite  large  compared  to  the  weaker 
intermolecular  cohesive  energies  associated  with  organic 


molecular  crystals  (which  are  due  mainly  to  van  der 
Waals  interactions),  since  it  is  a  result  of  is  numerous 
long-range  electrostatic  interactions  of  ionic  partners 
within  the  crystal.  This  value  can  be  directly  calculated 
by  simply  adding  up  all  interatomic  interactions  within 
the  crystal  lattice  if  the  microstructure  of  the  ionic  crystal 
and  a  reasonable  description  of  the  interatomic 
interactions  are  known.  Unfortunately,  these  two 
requirements  are  not  always  available;  further,  it  is 
preferable  to  have  a  methodology  for  predicting  this  value 
that  does  not  rely  on  experimental  data.  Toward  this  end, 
Jenkins  et  al.  (Jenkins  et  ah,  1999;  Jenkins  et  ah,  2002) 
developed  an  equation  that  correlates  the  inverse  cube 
root  of  the  molecular  volume  of  the  formula  unit  of  an 
ionic  crystal  MpX9  with  its  lattice  potential  energy 
Upot(M;,X?).  Jenkins  et  al.  parameterized  the  functional 
form  for  Upot  using  lattice  potential  energies  from  a 
previous  work  of  Jenkins  (Jenkins,  1998).  These  data 
consist  of  values  derived  from  thermochemical  cycles 
using  experimental  information  and  those  from  full  scale 
calculations  of  the  lattice  potential  energy  of  a  crystal  and 
the  information  used  in  the  parametrization  was  limited  to 
salts  containing  alkali  metal  and  alkaline  earth  cations  (Li, 
Na,  K,  Rb,  Cs,  Mg,  Ca,  Sr,  Ba).  Several  possible  sources 
of  errors  exist  for  this  approach  in  its  application  to  the 
high-nitrogen  compounds  we  are  evaluating,  including  the 
variable  reliability  of  fitting  data  (lattice  energies)  and 
assumption  of  the  ionic  volumes  used  in  the 
parameterization.  In  the  Jenkins  approach,  the  cations  are 
very  small,  single  spherical  elements,  whereas  many  of 
the  cations  in  the  high-nitrogen  salts  are  very  large. 

Because  there  is  no  way  in  which  the  Jenkins  et  al 
approach  can  be  assessed  directly,  we  can  only  compare 
the  summed  energies  of  the  Born-Haber  cycle  with  the 
few  experimental  data  we  have: 


Table  1.  Experimental  versus  Theoretical  for  Solid  Phase 
Heats  of  Formation. 


Species 

AHf,MX(s)  (kcal/mol) 

Experiment 

Theoretical3 

This 

Work 

l-5-Diamino-4- 

methyl- 

tetrazolium 

nitrate 

8.4 

41.7 

47.5 

l-5-Diamino-4- 

methyl- 

tetrazolium  azide 

138.1 

161.6 

165.4 

l-5-Diamino-4- 

methyl- 

tetrazolium 

dinitramide 

45.3 

92.1 

94.9 

a  Experimental  and  theoretical  column  values  from 
Galvez-Ruiz  et  al.,  2005. 
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Note  that  the  theoretical  result  generated  by  Galvez-Ruiz 
et  al.  (Galvez-Ruiz  et  ah,  2005)  invoked  a  Born-Haber 
cycle  using  reactions  along  with  proton  affinities,  whereas 
we  utilize  the  procedure  outlined  above.  Additionally,  the 
Galvez-Ruiz  et  al.  calculations  used  the  G2  method 
(Curtiss  et  ah,  1991),  which  has  different  empirical 
corrections  and  approaches  than  the  G3  methods  we  used. 
Note  that  the  two  theoretical  approaches,  which  would 
have  distinctly  different  sources  of  error,  are  very  close  in 
value  for  the  three  compounds.  They  both,  however, 
overshoot  the  experimental  values  by  30-40  kcal/mol. 
Since  the  common  element  between  the  two  theoretical 
methods  is  the  Jenkins  method  for  predicting  lattice 
energies,  we  suspect  the  source  of  error,  if  any,  is  in  this 
step. 

2.3.  Prediction  of  Aquatic  Toxicity 

There  exists  no  single  test  which  constitutes  a 
standard  indicator  of  aquatic  toxicity,  and  previous 
theoretical  and  experimental  studies  on  energetic 
materials  have  been  performed  on  a  variety  of  species. 
However,  the  Environmental  Protection  Agency  (EPA) 
has  greatly  facilitated  the  analysis  of  aquatic  toxicity  by 
developing  the  Fathead  Minnow  Acute  Toxicity  Database 
(EPAFHM).  This  freely  available  online  database 
encompasses  homogeneous-assay  Fathead  Minnow 
(Pimphales  promelas)  LC50  toxicity  data  collected  for 
over  600  industrial  organic  chemicals,  and  is  an 
invaluable  resource  for  the  QSAR  (quantitative  structure 
activity  relationship)  practitioner  (Russom,  2003).  We 
use  selected  data  from  the  EPAFHM  in  conjunction  with 
the  commercial  QSAR  package  CODESSA  2.7.2 
(CODESSA,  2005)  to  develop  a  5-parameter  predictive 
equation  for  Fathead  Minnow  LC50  toxicity. 

To  develop  this  expression,  87  compounds  were 
chosen  from  the  EPAFHM  according  to  a  few  key 
selection  criteria.  No  compounds  containing  halogens 
were  selected  on  the  basis  that  halogen-containing 
systems  present  a  different  mode  of  action  for  toxicity  and 
are  outside  the  realm  of  interest  for  these  systems. 
Compounds  containing  nitro  groups  in  as  many  different 
topologies  as  possible  were  chosen  to  capture  as  much  as 
possible  the  effect  of  nitro  group  placement.  Finally, 
emphasis  was  given  to  compounds  containing  nitrogen  in 
many  different  chemical  environments.  The  CODESSA 
package  was  used  to  perform  a  best  multilinear  regression 
analysis  on  hundreds  of  descriptors  (parameters)  to 
determine  the  best  expression  to  predict  Fathead  Minnow 
LC50.  Calculated  values  from  the  resulting  expression 
versus  experimental  Fathead  Minnow  LC50  values  for 
this  87  compound  training  set  are  shown  in  Fig.  3.  Results 
are  highly  satisfactory,  where  the  resulting  square  of  the 
regression  correlation  coefficient  R"  is  0.8407  (indicating 
a  high  degree  of  correlation  between  predicted  and 
experimental  values)  and  the  root  mean  squared  error 


(RMSE)  is  0.3901,  an  acceptable  margin  of  error  for  these 
purposes.  The  parameters  used  in  this  correlation  include 
the  octanol-water  partition  coefficient  (the  sole 
experimentally  derived  parameter  in  the  expression,  and  a 
common  component  of  biological  correlation  expressions) 
and  several  quantum-chemically  derived  parameters 
relating  to  reactivity  of  chemical  moieties  within  the 
molecule. 


Figure  3.  Calculated  vs  experimental  fathead  minnow 
LC50  toxicities  for  the  87-compound  training  set. 

It  has  been  assumed  for  the  purposes  of  this  project 
that  a  sufficient  level  of  accuracy  is  to  be  able  to  predict 
the  LC50  value  to  within  a  factor  of  three  (or 
equivalently,  a  LogLC50  value  to  within  0.5  log  units). 
The  results  shown  demonstrate  that  toxicities  for  the  vast 
majority  of  the  chosen  compounds  are  predicted  with  this 
level  of  accuracy.  Assessment  of  the  largest  outliers 
(compounds  whose  predicted  values  are  outside  of  this 
range  of  accuracy)  show  no  trends  in  chemical  group  or 
mode  of  toxicity.  This  QSAR  is  currently  being  tested  on 
known  and  notional  energetic  materials  to  validate  its 
inclusion  as  part  of  the  computational  chemistry  toolkit. 


CONCLUSIONS 

We  have  described  a  newly-developed  computational 
toolkit  composed  of  models  for  the  prediction  of 
properties  of  energetic  materials  related  to  performance 
and  environmental  hazard.  These  tools  can  be  used  in  the 
design,  characterization  and  assessment  of  a  material 
without  the  need  for  synthesis  and  measurement.  Our 
intention  is  to  augment  the  toolkit  with  a  variety  of  critical 
models  in  order  to  simulate  and  analyze  complex 
phenomena  across  multiple  time  and  length  scales. 
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